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Abstract. The magnetic correlations, local moments and the susceptibility 
in the correlated 2D Kondo lattice model at half filling are investigated. We 
calculate their systematic dependence on the control parameters JkI^ E'lid V jt. 
An unbiased and reliable exact diagonalization (ED) approach for ground state 
properties as well as the finite temperature Lanczos method (FTLM) for specific 
heat and the uniform susceptibility are employed for small tiles on the square 
lattice. They lead to two major results: Firstly we show that the screened 
local moment exhibits non-monotonic behavior as a function of V for weak 
Kondo coupling J/f . Secondly the temperature dependence of the susceptibility 
obtained from FTLM allows to extract the dependence of the characteristic Kondo 
temperature scale T* on the correlation strength U. A monotonic increase of T* 
for small U is found resolving the ambiguity from earlier investigations. In the 
large U limit the model is equivalent to the 2D Kondo necklace model with two 
types of localized spins. In this limit the numerical results can be compared 
to those of the analytical bond operator method in mean field treatment and 
excellent agreement for the total paramagnetic moment is found, supporting the 
reliability of both methods. 



PACS numbers: 75.20.Hr, 75.30.Mb, 71.27.+a 
1. Introduction 

The understanding of strong correlations in mixed valent and heavy f-electron 
compounds is mostly based on two generic models described by the Anderson lattice 
Hamiltonian or, in the special case of almost integer valence, by the Kondo lattice 
Hamiltonian [1] [2] . The latter model which will be the subject of this work represents 
an extreme limit where correlations in the f-orbitals are taken as infinitely large due 
to the Coulomb integral Uf t with t denoting the hopping energy of conduction (c) 
electrons. On the other hand the Coulomb interaction U of conduction electrons (and 
also the one between c- and f-electrons) is completely neglected in this model. The 
(7 = Kondo and Anderson models have the advantage of a simple and meaningful 
mean field solution [21 |3] with the constraint of only singly occupied f orbitals 
implemented by an auxiliary boson in mean field approximation. In the lattice this 
leads to hybridized quasiparticle bands with an exponentially reduced hybridization 
gap. Close to the half filled case with one c- and f- electron per site they have an 
effective mass m* ^ mi, much larger than the conduction electron band mass nib [I]- 
The mean field Kondo lattice model may be merged with local density band structure 
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calculations in the renormalized band theory |3] leading to a powerful method to 
calculate realistic Fermi surfaces for Ce and Yb heavy fermion compounds, including 
crystalline electric field effects. 

It has been known that a strong Coulomb repulsion among the conduction 
electrons significantly alters the electronic properties of a metal. First and foremost, 
the kinetic energy as measured by the (effective) band width will be reduced and 
eventually vanish at a metal-to- insulator transition. Second, the trend towards 
localization is accompanied by the appearance of magnetic correlations. 

The central focus of the present paper is the influence of the above-mentioned 
correlation effects on the screening of local moments. Of particular interest is the 
question how the Kondo energy scale is affected by conduction electron correlations. 
We adopt a two-dimensional Hubbard model at half-filling for the conduction electrons 
where correlation effects have been found to affect strongly the electronic properties. 
As the interacting conduction electron Hamiltonian, i. e., the two-dimensional 
Hubbard model cannot be solved exactly for an infinitely extended system we extract 
the relevant information from suitably chosen finite clusters. This approach seems 
justified considering the local character of the quantities to be investigated. The 
relevance of the cluster approach is assessed by comparing its results to the predictions 
of a constrained mean- field theory for an infinite system in the limit of extremely strong 
Coulomb repulsion. 

Attempts to include the effect of correlations between conduction electrons which 
may become important when the latter originate from d-orbitals have sofar been 
mostly limited to the Kondo impurity models using various analytical techniques 
like perturbation theory [5j |6] and 1/Nf expansion for small U, a Schrieffer- Wolff 
approach |7J and RVB ansatz [S] for large U, a scaling approach |S] and also numerical 
techniques like NRG [9]. In the impurity problem it was concluded that the Kondo 
energy scale T* increases with U [S] [6| . 

Concerning concentrated systems with magnetic moments at every lattice site the 
situation is rather controversial. In the case of non- interacting conduction electrons 
(U=0) the electronic properties will be determined by the competition between the 
energy gain due to (local) Kondo singlet formation and to developing (long-range) 
magnetic correlations. It is to be expected that the subtle balance between these two 
tendencies will be affected by conduction electron repulsion. 

The correlated Anderson lattice problem was treated with a Gutzwiller variational 
method [10] and the ID correlated Kondo lattice with DMRG approach [11]. It was 
further investigated with QMC simulations and a fermionic bond operator method [12] 
used before in the uncorrelated case [13]. The strongly correlated {U » t) 'Kondo 
necklace' limit of the Kondo lattice was treated using exact diagonalization [14]. It 
was suggested [10] that the Kondo scale actually decreases for increasing U, in contrast 
to the impurity results. The competition of singlet formation and induced inter-site 
RKKY coupling and its signature in thermodynamics was studied in Ref. jl5l using 
ED and FTLM methods for finite clusters of the 2D uncorrelated Kondo lattice and 
in Ref. [16] for ID chains in the strongly correlated limit. 

The present work which uses the unbiased numerical ED and FTLM methods as 
well as the analytical bond operator technique has two clear objectives: i) Firstly 
to investigate the characteristic T=0 on-site and inter-site correlations and local 
paramagnetic moment {ti'Xoc) systematically in the whole {Jk,U) parameter range 
of the correlated Kondo lattice model, where Jk is the local antiferromagnetic 
exchange coupling. In particular the total local paramagnetic moment is an excellent 
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measure for the formation or breakup of the Kondo singlet state as function of control 
parameters [TS] . We will show that it exhibits an unexpected non-monotonic behaviour 
that has not been reported before, ii) Secondly we calculate the finite temperature 
behaviour of susceptibility and specific heat with FTLM on a small tile of the square 
lattice. From the maximum position we may extract the characteristic 'Kondo' 
temperature scale T* and in particular we investigate its systematic dependence on the 
correlation strength U. We will resolve the ambiguity mentioned before and show that 
T* increases monotonically with U. The ED and FTLM calculations will be performed 
on small Kondo clusters with open boundary conditions which have a more realistic 
one-particle density of states than the one with periodic boundary conditions [15| . 

We stress from the outset that in the context of the FTLM calculations on a finite 
cluster the meaning of T* is that of an average energy of low energy singlet-triplet 
excitations obtained from the maximum in the temperature dependent susceptibility 
and specific heat. This is indeed the way how an estimate of T* in a real concentrated 
Kondo compound may be obtained. It is not the same as the genuine single ion Kondo 
temperature Tk in the continuum limit which is exponentially small compared to the 
hopping energy. This difference between lattice and single ion Kondo scales appears 
already within approximate analytical theories [2] • 

Nevertheless a comparison of ED numerical results for finite clusters with 
analytical results obtained by the bond operator method for the extended lattice and 
in the large U limit is useful and legitimate. Since the former are exact it allows one 
to judge how far the large U model approximation is justified and give an estimate 
of the reliability and accuracy of the bond operator approximation. This theory has 
previously been used for quantum critical properties of the localized model in the large 
U limit [T71 [TSl [H] . A generalized fermionic version has also been applied to the finite 
U case [13j . In this comparison we focus on the paramagnetic local moment size and 
on-site correlations as function of U and on the finite temperature susceptibility. 

We note that the investigation of finite Kondo atom clusters on nonmagnetic 
surfaces is also of potential experimental interest. Sofar single magnetic Co adatoms 
embedded in metallic Cu„ clusters ^D] and two site Co-Co Kondo clusters connected 
by Cu„ chains [5T] have been investigated. However in these cases of a metallic 
substrate it is the simple conduction electron density of the substrate surface that 
plays the dominant role for the Kondo properties. The present investigation would 
be relevant for correlated Kondo clusters weakly adsorbed on an inert insulating 
substrate without chemical bonding of substrate and Kondo atoms. This has, to 
our knowledge, not been experimentally realized until now. Nevertheless it is clear 
that the investigation of finite size Kondo clusters is of great interest in itself and not 
only in view of the bulk Kondo lattice materials. 

2. Model definition and single particle spectrum 

The correlated Kondo lattice model , sometimes called Kondo-Hubbard (KLU) model 
contains three terms. As in the non-interacting model there is the hopping term for 
conduction electrons Cia ('^ t) and the Kondo exchange term (~ Jk) to localized spins 
Si at every site. Their localization may be considered as result of their large repulsion 
Uf. As mentioned before the Coulomb repulsion U for conduction electrons is usually 
neglected in the lattice model, although well studied for the Kondo impurity case. It 
is, however, very interesting to include its effect, both because it is physically present, 
in particular in intermetallic 3d-4f compounds and also because it allows to study 
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Figure 1. Single particle DOS for the 2D n.n. tight binding model in the 
continuum limit (full line) with logarithmic singularity at Fermi level for half 
filling (e = 0). Number of states for eight site cluster is shown for OBC (full bars) 
and PBC (open bars). Left inset: Eight site cluster: Right inset : k- points for 
PBC. 

theoretically the continuous crossover from the metallic Kondo screening case to the 
insulating case of interacting local spin dimers. The model is given by 

^^KLU = X! "^L^JT + X! '^^T^-ii + Jk^t, ■ Si (1) 

where riia = cl^Cia- In 2D the first, n.n. hopping term leads to the single particle 
spectrum given by the dispersion ek — ~t{coskx + cos ky). In the continuum limit 
(A^ — >■ c») it leads to a well known density of states (DOS) function shown in Fig.[T]by 
the full line which has a logarithmic singularity a e = corresponding to half filling. 
As was done in Ref. [15] for the uncorrelated KL model (U = 0) we will use the exact 
diagonalisation with Lanczos method for finite square lattice tiles. For the systematic 
classification of finite size tiles on the square lattice we refer to [22] . Due to the large 
number of eight states per site {Sz — ±i} (g) {n = 0;?! = 1(tz = ±i);7i = 2} in the 
KLU model only tiles with eight sites can be used (see inset of Fig. [1]) . This precludes 
also finite size scaling such as may be performed in the ID case [23l[Tl]. Furthermore 
as in Ref. [H] we restrict to the half filled case with the number of conduction electrons 
given by ric/N = 1 (N=8). Due to small cluster size (N) ric can have only discrete 
values and therefore the next possible value ric/N = 3/4 in the singlet spin sector is 
already far from the interesting half filled case. 

As discussed before |15| one may choose both open and periodic boundary 
conditions for the tiles (OBC and PBC respectively). Their single particle spectrum 
differs greatly. While for PBC most of the states are put directly at the Fermi energy, 
giving little resemblance to the continuum DOS the spectrum is more realistic for 
OBC with more even distribution of energies (Fig. [T]). It was shown in the case 
U = [15] that this leads to a qualitatively different moment screening or singlet 
formation for small Jk- In this work we use exclusively OBC. Naturally, since there is 
a finite gap between occupied and unoccupied conduction electron states for half filling 
(Fig. [1]), the true continuum behavior of the infinite Kondo lattice cannot be literally 
described, the numerical ED results presented in the following rather describe finite 
'Kondo molecules'. As mentioned in the introduction they may in fact be artificially 
generated by adsorption of 3d clusters on nonmetallic surfaces. On the other hand we 
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Figure 2. (a) Effective paramagnetic total moment squared as function of U/t for 
various Jx/i- A fiat maximum of tlie moment appears for weak Kondo coupling, 
(b) Corresponding varation of {iJ-i^^) as function of JkI^ for different [/. 

will focus mostly on the local moment and on-site or n.n. correlations. At least in 
the strong correlation limit {U jt ^ 1) it turns out that the local quantities obtained 
numerically for the eight-site cluster agree rather well with the infinite lattice results 
(Sec. ED. 

3. Numerical determination of the local moment and thermodynamic 
properties 

To determine the ground state and thermodynamic properties of our system, 
we diagonalize the Hamiltonian iteratively on a small tile, applying the Lanczos 
algorithm. We also use the low-lying states and their energies to construct 
thermodynamic expectation values, in particular to calculate the partition function 
and its derivatives |24]. We stress that this unbiased method is exact and reliable 
for finite clusters and can therefore be used as a reference to compare with other 
numerical or analytical methods. We note that the large number of (eight) states 
per site prevents the use of larger clusters and the application of a finite size scaling 
procedure such as described in detail in Ref. [22] for a two-dimensional 5 = 1/2 model 
with just two states per site. Therefore using the Lanczos method for the present model 
is limited to tiles with eight sites (a maximum of 10 sites seems attainable). Even then 
the largest subspace with total spin 5**°* = has dimension 739 162 (38 165 260 for 
iV = 10 sites). Nevertheless, as the comparison with analytical results will show, the 
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local moment and on-site correlations for the small tile give a realistic approximation 
to the extended Kondo lattice behavior. 

The total paramagnetic local moment is given by 



It contains the spin fluctuations of both the itinerant and localized spins and their 
antiferromagnetic on-site correlations due to the Kondo coupling. This is a central 
quantity for investigating the influence of correlations on the Kondo effect because 
the Coulomb repulsion U tends to localize the t spins and thus influences strongly 
their singlet formation with S spins. For the uncorrelated case {U — 0) it was shown 
previously that the on-site Kondo coupling also induces effective RKKY inter-site 
coupling between the localized spins. This mechanism is still present for non-zero U 
but it has to compete with the effective superexchange. Therefore, in addition to the 
total moment we will also investigate the on-site singlet correlation 5x3 = {^i ■ Si) 
and nearest neighbor (i, j) correlation function Sij — {Si ■ Sj). Due to the smallness 
of the cluster it is only reasonable to calculate it up to next nearest neighbors. The 
same would be true for the correlation function {Si ■ Tj) characterizing the 'screening 
cloud' of a given localized spin Si which is not evaluated here because its characteristic 
length scale is much larger than the cluster size. In the case of a single Kondo impurity 
in the continuum limit it extends over a range of ^ = vp/Tk where Vi? is the Fermi 
velocity and Tk the single ion Kondo temperature. Note however that even in this 
simple case the screening cloud and its length scale ^ have not yet been observed in 
reality (25! . 

With the FTLM approach [53] it is possible to calculate the temperature 
dependence of various important thermodynamic quantities like specific heat C (T) = 
NAkBC{T) and uniform susceptibility x(T) = NAti'Qigi^BY'xiT) per mole as 
cumulants of simple operators, e.g., 



where A'a is the Avogadro number, /xq is the magnetic permeability, /ib and fee are 
the Bohr magneton and the Boltzmann constant, g is the gyromagnetic ratio, and 
A'^ is the size of our tiles (or number of sites). For simplicity, we assume the same 
gyromagnetic ratio for both the t and the 5' spins. The thermal expectation values 
are to be calculated with the statistical operator p{T) = exp(— /377klu)j /3 = 1/ (fceT). 
In this work we restrict to zero magnetic field, and as we are working with a finite 
system size N , we can safely set {Sl°^)rp = in the first equation. Numerically, for 
each set of Hamiltonian parameters {t,U, Jk}-, the calculation of thermal averages 
involves two averaging procedures: Firstly, a random starting wave function is 
chosen, and the Lanczos algorithm yields up to 0(100) extremal eigenvectors and 
eigenvalues. Secondly between 100 and 500 of these Lanczos runs with different 
starting wavefunctions are performed. All eigenvalues and eigenfunctions obtained 
in this way are subsequently used to calculate the traces over the statistical operator 
as indicated in Eqs. (O, leading to the desired thermal expectation values. For details 
we refer to Ref. [24] . 



(mL) = (Sf) + (rf ) +2{S.. n) . 



(2) 




N 




(3) 




Figure 3. (a) On-site Kondo-singlet formation between localized S and itinerant 
T spins as function of Kondo coupling for various strengths of Coulomb repulison. 
(b) Nearest neighbor AF correlations between localized spins as function of U/t 
for various Kondo coupling strengths. For U=0 the n.n. coupling is of induced 
RKKY type while for larger U it is dominated by the superexchange term. 



Although we do not consider the magnetization M (per site) of the correlated 
Kondo model explicitly, we may get some information about finite-field properties by 
calculating the third-order susceptibility defined through the expansion 



(3) 



for an applied magnetic field B. Here x'^^K'^) 
higher order cumulant according to 

1 1 



(4) 

(7VaMo)(.9Mb)^X^^H7^) is given by a 



SI 



(5) 



This quantity is a measure of the nonlinearity of magnetization at low field strength B. 
It has been discussed previously for a localized spin model [5^ and plays a significant 
role in the discussion of some heavy fermion compounds j27| . 



4. Discussion of numerical results 



First we discuss results for the total local paramagnetic moment presented in Fig.[21[a). 
Two counteracting trends determine its size: On one hand the increase in U localizes 
the T -spins and leads to an increase of {tD from 3/8 for U = to 3/4 for 
U/t 3> 1. On the other hand for finite Jk the Kondo term establishes the AF 
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on-site singlet correlation (xj - Si) < 0. Therefore in the case Jk—0 the moment 
increases monotonically from 9/8 = 1.125 to 2 • (3/4) = 1.5 while for any finite Jk in 
the large U limit when becomes localized the moment will decrease due to singlet 
formation. For moderate 3k there is an initial increase of {^J^f^^) with U due to the 
first correlation effect and eventually a decrease due to the effect of Jk when r spins 
become localized at larger U. In between a maximum in (/^foc) develops as function 
of U. For larger 3k the singlet formation effect is so strong that it overwhelms the 
increase in (rf) and therefore no maximum appears beyond JA'/t=0.1. This effect 
is completely dominated by local correlations and should not be strongly influenced 
by the tile size. The corresponding local moment dependence on Jk for various U is 
shown in Fig. [2][b) . For the uncorrelated U = case the singlet formation leads to a 
continuous decrease of (/ifoj,) with Jk- Increasing U facilitates this formation due to 
the localization of r spins and the decrease becomes progressively steeper as function 
of Jk- In the limit [/ — > oo an arbibtrary small Jk will lead to the singlet ground 
state. 

The singlet formation may also be monitored directly by the on-site correlation 
(ti • Si) which is shown in Fig. ISja). Starting from zero at Jk = it becomes 
increasingly antiferromagnetic for growing Jk until the singlet value —3/4 is reached. 
Again the latter is approached more rapidly with increasing correlation strength U- 
The nearest neighbor induced AF spin correlations are presented in Fig.[3l^b). For U=0 
these correlations are of the induced RKKY type [15] . They decrease with increasing 
Jk because their evolution is impeded by the increasing singlet formation. For a fixed 
but small Jk the antiferromagnetic inter-site correlation first increases for small U 
due to the reduction of single and double occupancies and for larger U it falls off again 
due to the reduction of the superexchange At^ /U with increasing U. For larger Jk this 
effect is more pronounced as function of U. 

A recurrent topic in previous analytical theories of the correlated Kondo impurity 
model is the U-dependence of the Kondo energy scale. From a practical viewpoint it is 
frequently taken as the maximum position T* of magnetic specific heat or susceptibility 
which corresponds to an average singlet-triplet excitation energy. In the present 
context of finite size tiles the U or Jk dependence of the characteristic temperature 
scale T* can be conveniently obtained from the FTLM results according to Eq. ([3]) for 
both quantities. Although it is not identical to the single ion Kondo temperature Tk 
of the impurity model one may expect a comparable qualitative dependence on U or 
Jk may be comparable. 

The contour plot of Cv{T) in the T-U plane is shown in Fig. |lUa),(b) for 
Jx/t = 0,2 respectively. The former corresponds to the case of the pure Hubbard 
model. It illustrates that the high temperature peak at J7 = due to uncorrelated 
charge fluctuations splits into a lower temperature spin fluctuation peak associated 
with J = At^ /U and a broad continuum at very high temperatures. The charge 
fluctuation peak at U=0 correspond to the excitation energy of ~ Q.lbt between 
the highest occupied and lowest unoccupied states around e = in Fig. [TJ The 
behavior as function of U may be viewed as a precursor of the thermodynamic Mott- 
Hubbard transition for a finite tile size. When Jk is turned on the specific heat peak 
is dominated by the lowest singlet triplet excitation energies. It increases linearly for 
small U (Fig.lU^b)) and reaches a plateau around U ~ Qt corresponding to the highest 
excitation energies in Fig. [T] 

To get rid of the influence of charge fluctuations we also calculated the spin 
susceptibility with FTLM. The results are presented in Figs. [5lja)i6] as function of 
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Figure 4. Specific heat from FTLM as function of temperature and interaction 
strength U. (a) for the eight site Hubbard tile (Jx = 0) the C(T) = {NAkB)C(T) 
evolution as function of U shows the splitting of charge and spin fluctuation peaks 
(upper and lower peaks respectively) for U > it K, W/2 ( ly = 8t is the total 
band width), (b) for finite Jk the U=0 peak defines the Kondo temperature scale 
due to singlet formation. It increases with U until it reaches a plateau for U > 4t. 



Jki U respectively. In the former a contour plot shows the evolution of susceptibility 
maximum or Kondo temperature scale T* (black dots) as function of Jk for constant 
U /t = A. For comparison the T* maxima for the uncorrelated U = Q Kondo lattice tile 
are also given in Fig.[5l[b) by the crosses. In both cases T* increases linearly with Jk in 
the strong coupling limit Jk ^ 1 . For small Jk the values for ?7 = are considerably 
below the results for finite U which indicates the increase of the Kondo temperature 
scale T* with U. This is also seen directly in Fig.[6]where T* (black dots) can be seen 
to increase linearly with U for small U and reach a plateau for U > Qt similar as for the 
specific heat peak before. In that figure we also included the peak position of the third 
order susceptibility x*^'^^ (T) (white dots) which also increases with U. This quantity is 
experimentally accessible [27]. It peaks at a systematically lower temperature than the 
first order susceptibility. The reason is that it characterizes the strongest nonlinear 
increase of the local moment with field which happens precisely in the temperature 
region where the screened local moment and susceptibility at zero field drops to zero 
(see Fig. [SJ. For comparison the size of the spin gap from analytical calculation in 
the large U limit (Eq. (fTSl) ) shown by white crosses is seen to be sandwiched between 
the above FTLM values for the first order (black dots) and third order (white dots) 
susceptibility peak positions. 

5. Bond operator treatment of the strongly correlated 'Kondo necklace' 
limit 

A deeper insight into the phases and excitations of low dimensional quantum magnets 
requires the application of both numerical and analytical techniques. As demonstrated 
e.g. by the spin wave analysis of the frustrated 2D spin systems ( Ref. [21] and 
references cited therein) analytical results for the extended system, even if approximate 
or only available in limiting cases , are very helpful to understand the systematics of 
numerical ED results for finite clusters. It is therefore perfectly legitimate to proceed 
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Figure 5. Susceptibility from FTLM as function of temperature and Kondo 
coupling Jk- (a) Contour plot of dimensionless tx{T, Jk)- Tlie maximum (black 
dots) defines tfie Kondo temperature scale T* of singlet formation for U/t = i. It 
increases monotonically with Jk- (b) For comparison T* for the noninteracting 
iU = 0) case is shown by white crosses demonstrating the increase of T* with U . 



in a similar way for Kondo lattice models. Because we focus on the paramagnetic 
phase we will, however, use the bond operator approach as the appropriate analytical 
technique for comparison in the large U limit. 

In the following we will derive such analytical results in the limit of large 
conduction electron correlations {U ^ 2zt, z=4 is the coordination number) at half 
filling {ric = !)• In this case the charge fluctuations may be eliminated from the 
hopping and Hubbard terms leading to a pure exchange term of the now also localized 
conduction electron spins r^. Thus the appropriate large U limit of the model in 
Eq. ([1]) is a localized 2D spin Hamiltonian [55] given by 

HkN = j'^Ti ■ Tj + Jx'^Ti ■ Si (6) 

(ij) * 

Here the kinetic hopping term in Eq. ([1]) is replaced by an effective inter-site spin 
exchange J = ^ in the strong correlation limit of conduction electrons. The above 
Hamiltonian is of the generalized 'Kondo-necklace' type (in 2D) originally studied for 
ID chains with only the xy inter-site terms included It has, however, later been 
extended to higher dimension and including all components in the intra- and inter- 
site exchange with possible uniaxial anisotropics for both [171 IHl US]- The model 
may also be viewed as an asymmetric bilayer Heisenberg model [30] with and Si 
spins residing on different layers and only the former coupled by inter-site exchange 
J. These are generic models to describe quantum phase transitions between a singlet 
('Kondo') phase favored by the second term and an antiferromagnetically ordered 
phase favored by the first term. The transition occurs when the control parameter 
J/ Jk is larger than a the value {J / Jk)c defining the quantum critical point (QCP). In 
2D we have in the present isotropic model {J/Jk)c = 0.88 [TBldn]. Such transitions 
are frequently found in f-electron compounds where the control parameter may be 
varied by pressure or doping (i.e., chemical pressure). The above model allows to 
study the characteristic quantum critical behavior around the QCP disregarding the 
charge fluctuations. It has been investigated using numerical methods like Monte 
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Carlo (MC) simulations [3TJ|35], exact diagonalization methods ^14^ and dynamical 
mean field theory (DMFT) and also analytical methods like bond operator approach 
in mean field [T71 [TSl [TH] or hard-core boson treatment [321 ■ I'^ this approach the 
Kondo necklace Hamiltonian in Eq. ^ is mapped to a model of interacting singlet (s) 
and triplet (iQ,a = x,y,z) bosons by the bond-operator transformation [34]. These 
bosons describe the singlet and triplet states |s) = s^\0) and \ta) = it|0) (a = x,y,z) 
of the pair of spins (r^, Si) coupled by Jk at every site i. In terms of bosonic operators 
the spins are given by 

'''i,a — 2i^^\^i,o: ~ ^i^a^i ~ "^^a/Sjilj^ti^j), (7) 

where a, /?, 7 = x, y, z and eap'y is the totally antisymmetric tensor. The conventional 
spin and boson commutation rules are fulfilled. The restriction to physical states (only 
one boson per site) is expressed by the constraint s|si -I- J2a^i a^i.oi = 1 which may 
be implemented either on the mean field level or by hard core boson technique. The 
former is chosen here. 



5.1. Ground state energy and triplon excitations 

In mean field approximation the bond operator transformation leads to a bilinear 
bosonic Hamiltonian 

k.Q 

+ Ei'^^'^k(tLi-ka+*k„t-ka) (8) 
k 

where fj, is the chemical potential introduced to ensure the constraint and s — (s) is 
the mean field singlet amplitude. Furthermore we define ika = J2i exp(ikR„)tiQ 
and 7k = (cosfc^; -I- cos ky). Here we restrict to the nonmagnetic phase where the 
triplet amplitude t = 0. In this case s is always close to one. This Hamiltonian may 
be diagonalised by the Bogoliubov transformation 

Oka = COsh((/)k)ik,a + Sinh((/ik)ilka' 

alka = sinh((?!)k)ik,a + cosh((/ik)i^k,a' (9) 

with 

tanh(20k) — 

«k 

/k = ^Js^lu (10) 

1 1 _2 

dk = + -^Jk + -j^Js 7k 

Diagonalization leads to the bosonic triplon Hamiltonian 

k a—x,y,z 



Wk = \/di - 4/2 (11) 
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Figure 6. Susceptibility in {T,U) plane for fixed Jk- The Kondo temperature 
scale (black dots) increases first linearly with U and then reaches a plateau 
qualitatively similar to the specific heat result in Fig. [4] For comparison we 
also show the temperature of the x^'^H'^) maximum (white dots) and the triplon 
gap At (Eq. {TS}) of the large-U Kondo necklace model (black crosses). 



where Wk is the threefold degenerate {a = x,y,z) triplon dispersion and Eq{^,s) the 
ground state energy. Minimization of the latter leads to selfconsistency equations for 
the chemical potential and singlet amplitude s given by 



3 X — A d\^ 5 _2 

k 



27V^^ ^ 2 



3 J dk - 2/k 3 

k 

The smooth dependence of /i, s on the control parameter Jk / J is shown in the inset 
of Fig. [71 It also extends continuously across the QCP into the antiferromagnetic 
region [TO] . 

The bond operator method also gives an explicit expression for the singlet-triplet 
gap which is defined as Af = with Q = (tt, tt) [18 . Using Eqs. (|10llll) and setting 
J = At^ /U we obtain 

This expression for the gap may be compared to the Kondo temperature scale T* from 
the susceptibility maximum obtained in FTLM finite cluster calculations. 



5.2. The paramagnetic effective local moment 

In the numerical calculation the central quantity in the correlated Kondo model is the 
local paramagnetic moment (/ifo^.) 'v^hich is a direct measure of the singlet formation 
as function of Jx/t and U /t. This also extends to the large U limit where the local 
moment may be calculated analytically within bond operator approach. To gain a 
better understanding of the numerical results a comparison with the analytical method 
for large U is helpful. In spirit this is similar to the comparison of ED results for the 
Ji — J2 Heisenberg model with analytical spin wave calculations [22] . The elementary 
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excitations here are, however, gapped triplon modes in the paramagnetic ('Kondo 
singlet ') regime rather than spin waves in the antiferromagnetic broken symmetry 
state. 

First we transform the expression for {n^^^) given in Eq. ^ to bond operator 
basis. Using the defining relations in Eq. ([7]) with 6^ — \{t^ x t)^ we derive the 
operator identity 

ML-4^0t0„; {^il,) = l2{^%) (14) 

a 

where we used the spin space isotropy or degeneracy of triplon modes for the 
expectation value. Using the explicit expressions for Oa, defining Wk = cosh(0k), 
i>k = sinh(0k) and performing a Fourier transformation we obtain a concise form of 
the moment: 

k k k 

(15) 

The positiveness of the moment squared is ensured by Schwartz' inequality. Using 
Eq. ([TU)) we can write 

ukt'k = — ; Uk + «k = — (16) 
Wk Wk 

Inserting these explicit expressions into Eq. ((T5|) we obtain the final result for the 
paramagnetic moment as 

/;,2 , , 3r 1 ^ dk + 2/k .i 1 ^. dk-2/k .i 
Kl^ioc) 2VN^^d^-2fJ N^^d^ + 2fJ 

(17) 

where dk, /k are defined in Eq. (|10p . From this closed expression the moment may be 
obtained by the 2D momentum integration. 

It is worthwhile to consider the expression of (/^foc) strong Kondo coupling 

limit J/Jk < 1 or {U'^/UJk < 1). We define ek = ^ and use Wk = rfk\/l - 4- 
Expanding Eq. ([2]) in terms of ek we get 

(18) 

k k 

Which may be further evaluated to give 

Where in the last expression we replaced J — This equation demonstrates that 
the local moment is decreasing with increasing Jk (note that /x and s also depend on 
Jk / J ) and with increasing correlation U. This is what the full numerical calculation 
of Eq. (fTT]) discussed below indeed confirms. 
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5.3. Spin correlations and high temperature susceptibility 

To obtain a more detailed insight into the Kondo singlet formation and the influence 
of Coulomb correlations on it we have previously also calculated the evolution of 
on-site and next-neighbor spin correlations in the ED approach. It is desirable to 
calculate them with bond operator approximation in large U limit for comparison 
with numerical results. First we consider the on-site spin correlation function 
SKs{i) — {Ti ■ Si) between localized and itinerant spins and also the complementary 
partial local moments (t^) and {Sf). Transforming to bond operator representation 
and using the isotropy we obtain 

SKs{^) = \ [(a^L) - 3s-2((i„ + 4)2)] (20) 

(T-f) = (Sf) = |[Woc) + 352((t„ +4)')] 

The equality — (Sf) is only valid in the localized limit but does not hold in the 
original KLU model defined in Eq. ([ij for small U. The site index i is suppressed i.f. 
because these local quantities are uniform. Expressing the triplet operators t^ , tl. in 
terms of triplon eigenmodes by using Eq. Q this leads to 

((i.+iin^^E^^k + i'k)' (21) 

k 

from Eq. we finally obtain 

s.s 4[<,L>-34i:(|^)'i 

(x?)H«?)4[(.L) + 34E(^)'1 P^) 

If we take the sum of these expression according to Eq. ([2|) one obtains again the total 
local moment. 

The total moment in the large-?/ (localized spin) limit calculated in Sec. 15.21 is a 
zero temperature quantity determined by quantum fluctuations in the ground state. 
On the other hand the effective moment of a localized spin system is obtained from the 
high temperature behavior of the susceptibility. The latter may be calculated from a 
high temperature expansion. In this section we want to investigate how these effective 
moments are related and how well the high temperature expansion agrees with the 
FTLM results in the large-U limit. First we give the high temperature expansion for 
the total local moment at finite temperature T: 

(/xf„,)(T) = (S^)^ + (rf)^ + 2 (5. • T.)^ (23) 

where 

{A)j, = 1 tr {Ae'^"^^) , Z = tr (e"''^™) 

is the thermal expectation value of A formed with the Kondo necklace Hamiltonian. 
Expanding the statistical operator for large ksT we obtain 

2' 



(/xf„,)(T) = 25(5 + 1) 



1 



5(5 + 1) Jk ( Jk y 
3 fesT \kBTj 



AL(T) + 0(/|) . (24) 
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Figure 7. Paramagnetic local moment and its square in the large U-limit as 
function of J/Jk where J = At^ /U . Thick line gives the result of bond operator 
theory. Thin and broken line are ED results for the eight site cluster in the large 
U limit. The good agreement suggest that finite size effects for the local moments 
are small. The inset shows the dependence of singlet amplitude s and chemical 
potential /i in mean field bond operator theory (Sec.O. 



The prefactor is the paramagnetic moment of uncorrelated spins, where the factor 2 is 
due to the presence of two spins r and 5' per site with t = S — \/2. The parenthesis 
gives the first order high temperature correction of the local moment which depends 
only on the local Kondo exchange Jk- Likewise we can calculate the high temperature 
uniform susceptibility, which was defined previously as 

where the same g-factor has been assumed for both r- and S-spins. Evaluating the 
thermal average in high temperature approximation as before and using t — S one 
obtains for the high temperature susceptibility 

1 1 



xiT) 



3kBT 
1 



(26) 



1 - 



5(5+1) Jk 5(5 + 1) zJ 



1 



1 



ksT 3 
5(5+1) zJ 
3 



Formally, we can regard the last expression as the first two terms of a high-temperature 
expansion of a Curie- Weiss type susceptibility 

Xo{T) 



X{T) = 



(27) 



1 + e/T 

with an effective noninteracting Curie susceptibility Xo{T) and a Weiss temperature 
Q given, respectively by 



xo(r)-^r^ xAL(r), 



fcee 



3 fceT 
5(5 + 1) 



zJ 



(28) 
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Figure 8. Comparison of partial moments (top) and on-site AF correlations 
^Ksi^) = (■'"i ■'S'i) (bottom) from ED (full lines) and bond operator theory (BOT) 
(dotted line). The results for Sks have been scaled by the effective coordination 
number of a given cluster site {z^ff = 1..4; for BOT Ze// = 4). The deviations in 
both quantities have different signs and almost cancel in the total local moment 

i.e. /cbO = J for 5* — — 4. We write the high-temperature behavior 

of the susceptibihty x{T) in this suggestive form to stress that the temperature- 
dependent effective moment screening is due to the Kondo coupling Jk , and the Weiss 
temperature Q due to the effective magnetic intersite exchange J. One may conjecture 
that this expression may be vahd beyond its formal expansion regime. This may be 
checked by comparing Eq. ([26|) with the unbiased FTLM results for x{T) and jj,f^^{T). 



6. Comparison with numerical results 

It is an interesting question whether the previous analytical mean field results for the 
large U limit can be qualitatively compared to the numerical ED results for finite tiles 
in Sec.m The most obvious quantity to check this is {iJ^i^^)- The comparison is shown 
in Fig. [T] The inset of this figure gives the dependence of fj, and s on J/ Jk as obtained 
from the self-consistent solution of Eq. ([T^ almost up to the quantum critical value 
{J/Jk)c = 0.88 where the singlet-triplet spin gap vanishes and AF order would set in. 
We stay in the paramagnetic parameter range throughout this work. The AF region 
of the 2D Kondo necklace region has been explored in Refs. ]TE[ 119) . 

The main Fig.[7]shows the dependence of the local moment (a*?oc)^ ^^'^ square 
on J I Jk from bond-operator mean field theory (full black line) and from numerical 
results for the eight site cluster model in the large U limit with J — At'^/U. Here 
large U means large as compared to the tight binding bandwidth W = 2zt = 8t. 
The numerical results for various U jt in this limit are given as dashed lines. The 
agreement of analytical and numerical results is surprisingly good. This proves two 
points: Firstly the mean-field bond operator technique gives reliable on-site properties 
in the ground state. This is in agreement with the observation [19^ that values for 
{J/ Jk)c that determine the quantum critical point between Kondo singlet and AF 
phase are well reproduced by that theory. Secondly the agreement points to the fact 
that the numerical finite size effects on a local quantity like {iJ^i^^) are apparently 
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Figure 9. Susceptibility and effective moment from FTLM (upper full line). 
Dashed line is obtained from high temperature expansion formula using fcsS = 
0. 03794 as a fit parameter. It is considerably smaller than in Eq. II26I I due 
to reduced coordination number in finite clusters. Lower full lines from high 
temperature expansion with given by Ea.l l28l l using the approximate expression 
for x{T) in Eq. l|26]l . 

quite moderate. Noticeable deviations in the numerical and analytical results appear 
at larger J / Jk when the quantum critical point to AF order is approached or when 
U/t becomes too small. 

The comparison may also be made for the on-site spin correlation functions SKs{i) 
and the partial moments (r^) and {S'^) which are equal in the large U limit. The 
latter is shown in Fig. [8^ with numerical results from ED giving the proper local 
moment value 3/4 and the results from bond operator theory lying below. Their 
difference increases with increasing J, i.e. decreasing U. On the other hand the on- 
site AF Kondo correlation presented in Fig. shows the analytical result lying above 
the numerical ED value by a similar amount. The latter depends on the effective 
coordination number of the considered site in the finite cluster. The total moment 
(/^foc) sum of these individual contributions and because of the opposite sign 

of the differences the add up to zero approximately. This explains why numerical ED 
and analytical bond operator results for (iJ^i^^) in Fig.[7]agree so well despite the small 
cluster size. Including the previous results 1181 119| on quantum critical properties of 
the Kondo necklace model we can conclude that the energetics and local correlations 
are well described by the mean field bond operator method. However, it should be 
noted that the description of inter-site correlations and their U dependence is well 
beyond this approximation. 

Finally we come to the high temperature expansion for local moment and 
susceptibility in the large U limit and its comparison to the FTLM results. According 
to Eq. (|26l) the ratio Tx{T) / iif^^{T) should be proportional to the Curie Weiss factor 
in this equation. This comparison is shown in Fig. [HI It is seen that the result of the 
high temperature expansion lies considerably below the curve obtained from FTLM. 
Part of this discrepancy may be due to the fact that in FTLM the effective (average) 
coordination number in the eight site tile is smaller than z = 4. Therefore, in order to 
compare with high temperature expansion we may consider an effective Curie Weiss 
temperature QcS as fit parameter in Eq. (pS)) . Then the temperature dependence of 
X(T) from FTLM is well reproduced by the form of the high temperature expansion 



(7 / r = 20, JkI t=i 
J /t = 4t/U = 0.2 
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results in Fig. [HI 

7. Summary and conclusion 

In this work we investigated the local moment screening and spin correlations and 
in particular thermodynamic properties of the correlated Kondo lattice or Kondo- 
Hubbard model. We have focused primarily on the systematic dependence of local 
moment, spin correlations and susceptibility on the control parameters Jk and U 
and secondly on the U-dependence of the characteristic temperature scale T* . We 
used unbiased and exact numerical techniques like ED and FTLM for small tiles for 
the whole range of U /t and Jk/^ as well as approximate analytical bond-operator 
method for the large U-limit of the extended Kondo necklace model which contains 
only localized spins. In this limit both methods give excellent agreement on the 
dependence of the Kondo-screened total paramagnetic local moment (Fig. [7]) on J7 or 
J = At'^/U. Our ED and FTLM investigation may also be of particular relevance for 
finite nano- clusters of Kondo atoms adsorbed on surfaces. 

A first central result obtained for smaller U and Jk is a clear non-monotonic 
behaviour of {^J-f^^) on correlation strength is observed in the ED results. This non- 
monotonic U-dependence is also observed in the on-site and inter-site spin correlations 
where the latter are of mixed RKKY and superexchange character. It is the result of 
competing effects of conduction electron localization by U and local singlet formation 
due to Jk- 

The second central conclusion concerns the dependence of the Kondo temperature 
scale T* on correlation strength U for which controversial results have been reported 
previously. Our FTLM susceptibility and the analytical results presented in Fig. [S] 
show that it increases monotonically with U and reaches a plateau in the large U 
limit. In this limit the Kondo scale corresponds to the spin-gap for triplon excitations 
at the AF wave vector Q = (tt, tt). 

Furthermore a high temperature expansion for the susceptibility of the Kondo 
necklace model leads to an effective Curie- Weiss type expression where the local 
susceptibility is modified only by the Kondo term and the Curie- Weiss temperature 
is only due to the superexchange term. The resulting temperature dependence is 
similar to that of FTLM results in the large U limit and if the effective Curie- Weiss 
temperature is used as a fit parameter a quantitative agreement over large temperature 
region is obtained. 
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